Random Trimer Tilings 
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We study tilings of the square lattice by linear trimers. For a cylinder of circumference m, 
we construct a conserved functional of the base of the tilings, and use this to block-diagonalize 
the transfer matrix. The number of blocks increases exponentially with m. The dimension of the 
ground-state block is shown to grow as (3/2 1 ' 3 ) m . We numerically diagonalize this block for m < 27, 
obtaining the estimate Soo = 0.158520±0. 000015 for the entropy per site in the thermodynamic limit. 
We present numerical evidence that the continuum limit of the model has conformal invariance. We 
measure several scaling dimensions, including those corresponding to defects of dimers and L-shaped 
trimers. The trimer tilings of a plane admits a two-dimensional height representation. Monte Carlo 
simulations of the height variables show that the height-height correlations grows logarithmically at 
large separation, and the orientation-orientation correlations decay as a power law. 

PACS numbers: 



I. INTRODUCTION 

The study of statistics of densely packed polymers has 
long been of interest to physicists. Onsager had argued 
that solutions of long rod-like molecules should show ori- 
entational order at high densities 0. Flory's approxi- 
mate analysis suggested that linear rod-like molecules on 
a two-dimensional lattice should also exhibit an orienta- 
tional order at high densities p}. However, for the case 
of dimers on a lattice — the only case that is analytically 
soluble — it is known that for all non-zero monomer den- 
sities, there is no long-range orientational order In 
the limit of zero monomer density, one gets power-law 
decay of correlations for (bipartite) square and hexago- 
nal lattices 0, || , but only short-ranged correlations on 
the triangular lattice |(| . Recently there have been stud- 
ies of the dimer problem on the cubic lattice Q and of 
interacting classical dimers on the square Q and cubic 
9] lattices. 

Monte-Carlo simulations of Baumgartner show that for 
semi-flexible lattice polymers close to full packing there 
is no long range order, no phase transition, and the cor- 
relation length is of the order of the size of the poly- 
mer |Toj . However, the exact solution of a single semi- 
flexible polymer that is fully packed on the square lattice 
(the so-called Flory model) exhibits a low-temperature 
phase of crystalline order and an infinite-order transition 
to a disordered, critical high-temperature phase in which 
the critical exponents vary continuously with tempera- 
ture mi. 

It is generally believed, but not proved, that in the 
continuum case in three dimensions, long needle-like 
molecules would undergo a isotropic-nematic transition. 



* anandaOtheory. tifr .res . in 
t ddhar@theory.tifr . res. in 

* jesper .jacobsen@u-psud.fr 



In two dimensions, a spontaneous breaking of continu- 
ous rotational symmetry is not allowed, but there is a 
Kosterlitz-Thouless phase with power-law decay of orien- 
tational correlation functions |l2t 1 1 3| . But the situation 
is less clear for systems of hard-core molecules on a lat- 
tice. In the limit of high density, one can get a solid-like 
phase where one of the sublattices is preferentially occu- 
pied, e.g., in the cases of hard squares and hard hexagons 
hj. The same behavior is seen in Monte Carlo studies 
with lattice models of extended hard discs ^5|. However, 
triangular trimers on the triangular lattice can be solved 
exactly, and do not show a long-range order even at close 
packing 0. De Gennes has argued that long straight 
needles may not show an ordered phase on the square 
lattice There are not many studies of other molec- 
ular shapes in lattice models, as realistic modeling of ac- 
tual experimental system assemblies of different shaped 
molecules (e.g., ellipsoids, banana-shaped molecules, etc.) 
is better done in the continuum space. Tilings by re- 
shaped trimers and T-shaped tetramers of to x oo strips, 
for to < 5 have been studied using the transfer matrix 
technique earlier [TH Il9j| . 

In this paper we study random trimer tiling of the 
square lattice with horizontal and vertical trimers (Fig. 
0]). We wish in particular to assess if there is a long- 
ranged correlation of the orientational order in this prob- 
lem in the limit of zero monomer density. The problem is 
first addressed in the geometry of semi-infinite cylinders 
of size to x oo . We show how to set up the corresponding 
transfer matrix and numerically diagonalize it for to < 27. 
We use these data to extrapolate to the to — > oo limit, de- 
termining in particular the entropy per site in a random 
trimer tiling in the thermodynamic limit. We also prove 
the existence of a family of matrices which commute with 
the transfer matrix and show that the transfer matrix de- 
composes into a number of blocks which is exponentially 
large in to. Our results for the free energy and various 
correlation functions are consistent with a conformally 
invariant system of central charge c = 2.15±0.2. In 
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FIG. 1: A trimer tiling with periodic boundary condition in 
vertical direction. 




FIG. 2: Choice of unit vectors 

particular, the correlation length £ m is found to increase 
linearly with m. 

It is known that trimer tilings admit a two-component 
height field representation 20]. In the plane, we use 
this to study correlations of the height field by Monte- 
Carlo simulations, and find that the correlation function 
(h(r)h(O))) varies as log ( | r | ) for large |r|. We have also 
studied the orientational correlations in the simulations, 
and find that they decay as power laws. 

The plan of the paper is as follows: in Section [Q] we 
present the construction of a height representation for a 
trimer tiling; in Section II I II the constants of motion are 
obtained for trimer tilings on a cylinder; in Section llVl thc 
number of disjoint sectors is calculated by the generating 
function formalism; in Section E] we set up the transfer 
matrix for the problem, and in Section IVII we present 
the results of numerical diagonalization of the transfer 
matrix; in Section I Vlll we discuss the trimer-trimer cor- 
relation functions and in Section lYlIII height- height cor- 
relation functions, both obtained by Monte-Carlo simu- 
lations. Our conclusions are presented in Section llXl An 
appendix adapts the working of Section IIVI to a more 
general tiling problem. 



II. TRIMER TILINGS AND HEIGHT 
REPRESENTATION 

We represent each trimer as a set of three consecutive 
squares along a line, oriented horizontally or vertically. 
For each tiling of a two-dimensional plane by trimers (Fig. 
1), one can define a configuration of a height model at the 
vertices of the square lattice, where the heights are two- 
dimensional vectors 0,112 as follows: 



Choose two-dimensional vectors ei , e 2 and e3 such that 
they satisfy the condition ei + e2 + e3 = 0. A convenient 
choice is ei = (1,0), e 2 = |(— l,v / 3),e 3 = |(— 1,— \/3), 
as shown in Fig. [21 Equivalently the vectors can be repre- 
sented as complex numbers ei = 1, e 2 = lu and e3 = u> 2 , 
where u = e 27 ™/ 3 . The height h(i,j) at any site is 
an integer linear combination of basis vectors ei, e 2 and 

The lattice edges are assumed to be oriented rightwards 
or upwards, and they are labeled with vectors ei, e 2 and 
e3 periodically as shown in Fig. [31(a) and (b). The la- 
beling is such that moving along a horizontal or vertical 
line from any vertex up or right encounters a periodic 
sequence of labels ei, e 2 , e 3 . This rule still leaves some 
freedom in choosing the sequence of bonds; two conve- 
nient choices are shown in Figs. [3(a) and (b). 

Now, for any given tiling of the square lattice, the 
height field h(i,j) is defined so that if the directed edge 
£ from site a to the site b does not belong to the interior 
of a tile (i.e., it forms part of the boundary shared by 
two tiles), and it has label e Q , then h(h) — h(a) = e Q . 
This determines the heights at all vertices up to an unim- 
portant additive constant. The constant can be fixed by 
arbitrarily choosing h(0) = 0. 

Note that for the two choices of edge-weights shown in 
Fig- El the height difference along any edge of trimer has 
modulus 1. For (a), it value is v7 or 2 for the internal 
edges of trimers, and for the choice (b), it only takes the 
value v7. An example of the values of the height field 
following the convention (a) is illustrated in Fig. |^c) for 
a particular configuration of trimers. 

The choice (a) has a particular advantage. It can be 
shown that any trimer covering of a plane can be obtained 
from any other by a sequence of the basic flip operation, 
in which three adjacent horizontal trimers are replaced by 
three adjacent vertical ones (Fig. |3J- Under a basic flip 
move, it is easy to check that the height changes only at 
four sites and the modulus of the change in height || A(h)| 
is always 3. We can think of the values of the height field 
as forming a triangular lattice on the complex plane. This 
lattice can be broken into 9 sublattices (Fig. [5J , such that 
even after any such flip, the value of height stays on the 
same sublattice. Also, different sublattices of the height 
field are in one-to-one correspondence with a 9 sublattice 
decomposition of the original square lattice (Fig. [SJ . 



III. CONSTANTS OF MOTION FOR TRIMER 
TILINGS ON A CYLINDER 

In this section we consider the geometry where the 
square lattice has been wrapped on a cylinder of circum- 
ference m. The fully-packing constraint on the trimer 
tilings then implies strong constraints on the different 
configurations of trimers along a row of length m. In 
fact, given some local configuration along a row, many 
local configurations along a different row are disallowed. 
If we think of the row-to-row transfer matrix as an trimer 
evolution operator for configurations on a line these con- 
straints can be described in terms of some constants of 
motion under this evolution. 
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FIG. 3: (a) A trimer tiling of 3 x 5 lattice with assignment of 
labels to the edges of the lattice (b) An alternative assignment 
of labels (c) the height configuration corresponding to label in 
(b), where (m, n?) denotes the height h(ni, 712) = niei+ri2e2. 



FIG. 4: The flip move exchanging three vertical trimers with 
three horizontal trimers. 



The simplest constant of motion of this type can be 
constructed in terms of the invariants for loops p^ . 
Given a trimer tiling, define an allowed loop as a sequence 
of nearest neighbor bonds on the lattice that returns to 
the starting point, and does not intersect itself, and none 
of the steps crosses a trimer. In other words, an allowed 
loop goes along the boundaries of trimers. Define Q as 
the group generated by two generators a and b, which 




FIG. 5: (a) 3 x 3 superlattice decomposition of the square 
lattice into 9 sublattices. (b) Height variables corresponding 
to a tiling of the square lattice form a triangular lattice. 

satisfy 

a 3 b = ba 3 

ab 3 = b 3 a. (1) 

Now we attach weights (elements from the group Q) to 
each step of the loop. These read respectively a, a -1 , 
b, and b 1 for a step of the loop to the right, to the 
left, up, and down. Finally, to each loop L we assign a 
weight wl € G as the ordered product of the weights at- 
tached to the steps along the loop. Eqs. Q imply that 
wl does not depend on the way (starting point and direc- 
tion) in which the loop L is traversed when building up 
the product weight. Moreover, it is easy to see that for 
a non- winding loop (i.e., an allowed loop on the cylinder 
which is nomotopic to a point), wl is equal to identity. 
(The proof goes by induction on the number of trimers 
enclosed by the loop.) 

Also, Wl takes the same value for any allowed loop that 
winds around the cylinder. Indeed, let L be such a loop. 
Then it is easy to see that w l does not change if the loop 
is deformed locally so that the number of trimers below 
it changes by ±1. Thus, wl is a constant for the tiling. 
However, given two products wl and Wl> of the genera- 
tors a* 1 , b* 1 , checking whether wl — wl 1 by using the 
rules is nontrivial. We now describe a different con- 
struction that is equivalent to this, but more convenient 
to use. 

Fig.[S]shows a partial tiling of the plane, starting from a 
base (shown hatched in the figure) . Each trimer occupies 
three horizontally or vertically consecutive squares. Peri- 
odic boundary conditions in the horizontal (x) direction 



4 



are assumed. The base is supposed without overhangs, 
and can therefore be specified by its height profile, H (x). 
The only constraints about the partial tiling we assume 
is that it has no holes, and no overhangs. In particular, 
the tiled region can also be described by a height profile 
H{x). 

We now construct an invariant of the tiling as follows. 
First, let A, B and C be any three non-commuting ma- 
trices satisfying the condition 



bbbbbbbbbbbbbbb 



B 3 



Clearly such matrices exist; one family of possible choices 
is given by 



1 Ai ,ui 
u> v\ 
to 2 



B 



1 A 2 M2 

LO 2 V2 

to 
1 

A3 LU 

A*3 ^3 



(3) 



where u> = e 2 ™ 1 / 3 , and Ai, i>i are any complex numbers. 

Next, we assign to each square one of the three colors a, 
b or c, where each color corresponds to a value of its verti- 
cal coordinate y (mod 3) (see Fig. |HJ| . The height profile 
H(x) can then be characterized by a word l\l2 ■ ■ ■ lm over 
the letters a, b, c, where for any x = 1, 2, . . . , m the letter 
l x specifies the color corresponding to y — H{x). For ex- 
ample, the word characterizing the partial tiling in F'ig. [H] 
reads bbbaaabcccbccba. Correspondingly to this word, we 
construct a functional of H{x), denoted by J(H(x)), 
which is defined as the trace of a product of matrices 
A, B, C, where each factor is obtained from a letter in 
the word by replacing ai-^A,6i-^B,ci-^C. For exam- 
ple, for the sequence given above, we get the functional 



J{H{x)) = Tr [BBBAAABCCCBCCBA] 



(4) 



Note that by the usual cyclic properties of the trace, the 
functional depends only on the height profile, and not on 
the starting point. We can therefore write J(H(x)) = 
J(H). 

We now have the remarkable theorem: 



J(H) = J(H 



(5) 



In other words, J is the same for all valid partial tilings, 
grown from the same base, and is equal to its value for 
the base. 

The proof is by induction on the number of trimers in 
the tiling. It is clearly true for no trimers. If we add a 
vertical trimer at x, the height increases by 3, and the 
letter l x does not change. If we add a horizontal trimer, 
it must be done on a locally flat substrate. Then in the 
word, the substring aaa may be replaced by bbb, or bbb by 
ccc, or ccc by aaa. But by Eq. @, this does not change 
J. Q.E.D. 
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FIG. 6: Partial filling on top of a given base (shown hatched). 



For the particular choice of A, B, C given in Eqs. JSJ, 
we actually have the additional property A 3 = B 3 = 
C 3 = I. For any choice of matrices having this additional 
property the word UU ■ ■ - l m can be transformed into an 
"irreducible word" |2J| by recursively deleting any subse- 
quent (modulo m) occurrences of three equal letters [i.e., 
lili + iU + 2 = aaa, or bbb, or ccc] in the original word. In 
the example shown in Fig. El the base is characterized by 
the (reducible) word cccbbaaaccaaaba, which then corre- 
sponds to the irreducible word bbecba. By the theorem, 
this is also the irreducible word of any partial tiling built 
on that base. 

The invariant J is a polynomial function of the pa- 
rameters Ai , jii , Vi in the matrices © . We can therefore 
expand J in multivariable power series in these variables. 
By the theorem, all the coefficients in this expansion are 
the same for any allowed tiling, and are constants of mo- 
tion for the transfer matrix. 

For trimers tilings on a torus, we can define other in- 
variants J(n x ,n ) corresponding to other homotopy classes 
of non-contractible loops. The allowed homotopy classes 
are characterized by the winding numbers n x , n y along 
the two coordinate directions, with n x A n y = 1 and 
n x + n y > 1 (2£j. Our previous invariant is J = Jn,o)' 
However, the invariants 3(n x ,n y ) are n ot all independent, 
as can easily be inferred from the simple example when 
l7(i,o) corresponds to abababab . . .; and there is then only 
one possible trimer tiling which completely fixes the val- 
ues of all other loop variables. 



IV. DECOMPOSITION OF PHASE SPACE INTO 
DISJOINT SECTORS 

The word associated with the base can be any sequence 
of m letters chosen from a, b, c. The number of all pos- 
sible base profiles Ho (x) (modulo the addition of vertical 
trimers) for a cylinder of width m is therefore 3 m . 

We have already seen that two height profiles Hi(x) 
and H2(x) are reachable from each other if and only if 
J(Hi(x)) = J{H2{x)). Hence the transfer matrix for 
this problem has a block diagonal structure, with no tran- 
sition possible between configurations with different J . 

It is straightforward to determine the number of irre- 
ducible words of length n > 2 that start with any two 
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specified letters. Let us denote by N d (n) (resp. N s (n)) 
the number of irreducible words in which the initial two 
letters take fixed different (resp. same) values. Using the 
condition that in an irreducible word we cannot have any 
three consecutive letters being identical, it is easy to see 
that they satisfy the recursion relation 



N d (n + 1) = 2N d (n)+N s (n) 
N s (n+1) = 2N d (n) 



(6) 



for n > 2, with the initial condition N d (2) = N s (2) = 1. 
These equations are easily solved giving 



N d (n) = -1_[-(1-V3)» + (1+V3)»] 



(7) 



Ns{n) = ^-[(3 + V3)(l - V3) n + (3 - V3){1 + V3) n } 

The total number of irreducible strings t n of length n is 
then obtained by taking into account the multiplicities 
due to the possible choices of the two initial letters: 



t n = 6N d {n) + 3N s (n) . 



(8) 



The total number of sectors for a cylinder of width 

to (with 3|m) is then ^fe^o^ 3 ^ Note that this number 
is exponentially large in to. For comparison, for dimcr 
tilings (with 2\m) the number of sectors is just 2m + 1. 

Note that each of the t m sectors corresponding to an ir- 
reducible word of length to contains only one state of the 
transfer matrix (these sectors are "stuck"). On the other 
hand, the to = 1 sector corresponding to the empty irre- 
ducible word comprises a number of states in the transfer 
matrix which grows exponentially with to. 

It is straightforward to obtain precisely this latter num- 
ber, i.e., the number of words of n letters that reduce to 
the empty irreducible word. Namely, it corresponds to 
the dimension of the transfer matrix in the ground state 
sector, i.e., the sector from which the free energy is ob- 
tained. 

To this end, let us first define 

G a = aaa + abbbaa + aabbba + acccaa 

+ aaccca + abbbeccaa + accebbbaa + ... (9) 

as the formal sum over all unfactorizable words with ini- 
tial letter a that are reducible to the empty word. Here 
unfactorizable means that the words contributing to G a 
must not be the concatenation of two non-empty words 
each of which is in turn reducible to the empty word. 
Note also that the requirement that the initial letter be 
a implies, by the property of unfactorizability, that the 
last reduction before reaching the empty word is of the 
type aaa i— ► 0. We similarly define G b and G c . The sum 
in Eq. @ can then be expressed as 



G„ 



a[l + {G b + G c ) + (G b + G c ) 2 + ■■■] 
a[l + {G b + G c ) + (G b + G c ) 2 + ■■> 
1 1 



l-(G b + G c ) l-(G b + G c ) 



(10) 



Now substituting a = b = c = x in Eq. we obtain the 
generating function for irreducible words with a formal 



weight x per letter: 
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(x) = G a (a = x) g 3n i 



.3 it 



(ii) 



where g 3n is the number of different unfactorizable words 
that start with a given fixed letter and are reducible to 
the empty word. By Eq. I)10[l. g(x) then satisfies the 
equation 



l(x){l - 2g(x)} 2 =x 3 . 



(12) 



This is a cubic equation in g{x) and can be solved ex- 
plicitly. Among the three solutions for g(x), two can be 
discarded as unphysical on the ground that g(0) ^ 0. The 
last, physical solution can be expanded into a polynomial 



scries in x, as 



4a; 6 + 28a; 9 + 240a; 
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2288a; 



15 



23296x 18 + 248064a; 21 + 2728704a; 24 



Apart from the trivial root in x = 0, g{x) has two non- 



trivial coincident roots for 
x c , g{x) varies as 



„3 _ „3 



2/27. For x near 



g(x = x c -5) = l- A5 1 ' 2 + oiS 1 ' 2 ) . (13) 
6 



This implies that for large n 

'27\" I 



.9.3 r, 



.4 



,3/2' 



(14) 



Using g(x), we can now construct the generating func- 
tion H(x) of all words (factorizable or not, and with any 
initial letter) that are reducible to the empty word. We 
have clearly 

H = l + (G a + G b + G c ) + (G a + G b + G c ) 2 + --- 

= l-(G a + G b + G c ) (15) 
Setting a = b = c = x in Eq. fTKf as before, we get 



H(x) 



H 3n X 



3/ 1 



(16) 



where H 3n is the total number of words of length n that 
are reducible to the empty word. The leading terms of 
the polynomial series are 

H(x) = 1 (17) 
1 - 3g[x) 

= 1 + 3a; 3 + 21a; 6 + 183a; 9 + 1773a; 12 + 18303a; 15 
+ 197157a; 18 + 2189799a; 21 + 24891741a; 24 + • • • . 

For x near x c , H(x) w H(x c ) — A(x c — a;) 1 / 2 , whence H 3n 
also varies as 



H. 



-'!// 



A 



27 



1 



,3/2' 



(18) 



for large n, where A is the same constant as in Eq. I|14|) . 

The coefficients appearing in Eq. (|17(l coincide with the 
observed dimension of the transfer matrix in the ground- 
state sector (see Table [I] below). 

Finally, let us note that the working of this section can 
be adapted to a more general tiling problem. This is 
relegated to Appendix A. 
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FIG. 7: Labeling of the R matrix. 

V. SETTING UP THE TRANSFER MATRIX 

The number of ways N n (to) to tile a cylinder of width 
to (with 3|m) and height n, with periodic boundary con- 
ditions in the m-direction and free boundary conditions 
in the n-direction, can be found as (0|(T m ) n |0). Here, T m 
is the transfer matrix that adds one row of width to, |0) 
is an initial state corresponding to an initial horizontal 
base, and (0| is a projection operator on the state corre- 
sponding to a final horizontal base. For finite n, diago- 
nalizing T by means of a similarity transformation leads 
to an expression of the form N n (m) = J^. ai(m)[\i(m)] n , 
where Xi(m) are the eigenvalues of T. The correspond- 
ing amplitudes oti (m) depend on our choice of boundary 
conditions in the n-direction. We are mainly interested 
in the limit n>l, for which one has simply 

N n (m) « [Ai(m)] n (19) 

where Ai(m) is the largest eigenvalue of T m . The corre- 
sponding entropy per site is then 

S m = -\og\ x {m). (20) 

TO 

We shall now describe two different ways of construct- 
ing the transfer matrix. 

A. First construction 

The first construction is conveniently described in 
terms of a more general tiling problem (cf. Appendix A) 
of tiling the plane by horizontal p-mers (of size p x 1 ele- 
mentary squares) and vertical g-mers (of size 1 x q). The 
trimer case is recovered for p = q = 3. 

We first shift the tiling by half a lattice spacing, both 
horizontally and vertically, with respect to the underlying 
square lattice. The tile boundaries then intersect some of 
the lattice edges, and it is natural to describe the tiling 
by assigning an appropriate variable ('spin') s» to each 
lattice edge i. 

An edge i intersecting a tile boundary has Sj = 0. Each 
vertical g-mer encloses q — 1 lattice edges which are not 
intersected by its boundaries; the lowest of these edges 
has Si = 1, the second-lowest s,- = 2, and so on, and 
the highest edge has Sj = q — 1. Each horizontal p-mei 
encloses p — 1 lattice edges which are not intersected by 
its boundaries; the leftmost of these edges has s, = q, 
the second- leftmost has Sj = q + 1, and so on, and the 
rightmost edge has s$ = q + p — 2. 

Using a construction well-known from the theory of 
integrable systems, the row-to-row transfer matrix can 
be written 

T m = TTq Rom ■ ■ ' R02R0I > (21) 



where each of the matrices i?ofc act on two ingoing edges 
(labeled and k), joins them by adding one vertex, and 
produces two outgoing edges (labeled k' and 0'), as shown 
in Fig. [7| The trace over the edge corresponds, graphi- 
cally, to adding the first horizontal edge of a new row, and 
making sure that it joins to a horizontal edge carrying the 
same spin so, once the row has been completed. 

It remains to specify the elements of the matrix 
Rok(so, Sk', s' k , s' ). These are one for the cases 

r (0,0; 1,0), (0,0; q,0) 
(0, + 0) with 1 < i < q - 2 
(s ,s k ;s' k ,s' )= { (0,9-1:0,0) 

(i, 0; 0, i + 1) with q < i < q + p - 3 

{ (q+p-2, 0;0,0) 

(22) 

and zero for all other cases. Indeed, the first line in 
Eq. (|22|l corresponds to the lower left corner of any tile; 
the second line to the interior of a vertical tile; the third 
line to the upper right corner of a vertical tile; the fourth 
line to the interior of a horizontal tile; and the fifth line 
to the upper right corner of a horizontal tile. Finally, the 
trace in Eq. (|21() is over sq = 0, q, q + 1, . . . , q + p — 2. 

In some applications it might be of interest to give 
different weights to horizontal and vertical tiles. This 
can readily be done, by attributing the desired weight to 
the lower left corner of each tile, corresponding to the 
first line of Eq. 

The transfer matrix T m is constructed in the base of 
spins states (si, S2, . . . , s m ) corresponding to an initial 
row of vertical edges, with each Sfe G {Q, 1, . . . , q—1}. The 
factorization H21|) of T m is particularly suited for using 
sparse-matrix and hashing techniques, so that T m can 
be multiplied onto an arbitrary initial vector in time ~ 
mdim(T m ). The first few eigenvalues of T can be found 
by an iterative scheme (the so-called power method p6|) 
based on iterating such multiplications. 

The ground-state sector, corresponding to growing the 
tiling from a horizontal base, is obtained by choosing the 
initial vector so that the reference state (s\, . . . , s m ) = 
(0, . . . , 0) carries weight one, and all other states carry 
weight zero. The whole state space for the given sec- 
tor is constructed automatically in the iterative process. 
For the ground state sector, the total number of states 
is found to be given by Eq. I|17f) : this constitutes a use- 
ful check of the numerical algorithm. Excited sectors, 
corresponding to non-horizontal bases, can be similarly 
accessed by choosing another appropriate reference state 
as the initial vector. 

The above construction has enabled us to numerically 
diagonalize T m for m < 27. The corresponding eigenval- 
ues are tabulated in Table [I] 

B. Second construction 

The second construction of the transfer matrix is de- 
scribed here for the original trimer tiling problem. We 
first define the tiling at level n as the set of all tiles that 
have at least one square with y-coordinate less than or 
equal to n. And H n (x) is the height profile for this set 
of tiles above y = n. Then H n (x) clearly lies between 
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and 2, and may be characterized by a sequence of the 
type 21222000122110 • • • . The transfer matrix T c ,c is 1, 
if the height configuration C can be reached from C by 
adding some tiles, otherwise 0. 

The transfer matrix is constructed simply as follows: 
Let the configuration C be given by some sequence of 0's, 
l's and 2's as above. We add trimers to this configura- 
tion to all sites with height 0. We can place a horizontal 
trimer at any place where three consecutive sites have 
height 0. This changes the heights at these sites to 1. We 
place zero or more horizontal trimers this way. To each 
remaining site with height we add a vertical trimer, so 
that the height of that site becomes 3. Now there are no 
sites with height 0, and the maximum height is 3. Fi- 
nally we measure heights from a new reference point one 
unit higher, and decrease all heights by 1. This gives the 
new height configuration C with heights given again by 
a sequence of 0's, l's and 2's. 

There is a simple way to represent this transfer matrix 
as a spin hamiltonian. Denote the three possible heights 
at site i by a quantum spin that can be in any one of three 
orthonormal states |0), |1) and |2). We denote the three 
states at site i by |0j), |L) and |2,). Define the operator 
Si by 



Sflhi) = IM>> where h-ti = 1 (mod 3). 



(23) 



And we define Pf as the projection operator for the state 
|0j>, i.e. 



TO 



10,), P?\h) 



0. 



(24) 



Then it is easy to see that the transfer matrix can be 
written as 



T = Trn 



Sr Pf 
OOP? 
P? 



(25) 



To prove this, we only need note that expanding 
the product, the only nonzero terms are of the form 
S-S-P°P P S- . . ., where we have a string of S~ at 
consecutive sites, interspersed with the product of three 
P°'s at consecutive sites. 

As an example, let us consider random horizontal and 
vertical trimer tilings of a square lattice with periodic 
boundary conditions in the horizontal direction, i.e. , an 
infinite cylinder of width to. For convenience we choose 
to to be a multiple of 3. A typical tiling of the cylinder is 
shown in Fig. ^ It is easy to see that on a 3 x oo cylin- 
der a horizontal trimer can be followed by three vertical 
trimers, or by another horizontal trimer in exactly three 
ways due to the periodic conditions in the horizontal di- 
rection. The possible height configurations are thus '000', 
'111', and '222', and the to = 3 transfer matrix is given 
by 



We can reduce the size of the transfer matrix by work- 
ing in a sector where the basis vectors are invariant under 
translations and reflections. Thus for 6 x oo cylinder, in 
the sector where the irreducible word is empty, we have 
six basis vectors, '000000', '000111', '111111', '000222', 
'111222', '222222', and the other vectors are related to 
these by symmetry. The number of vectors needed for 
the transfer matrix for a 3n x oo cylinder is less than 
H^n by approximately a factor of 6n. The resulting size 
of the transfer matrix for n < 8 is shown in Table [I] 
By using the rotational and translational symmetries the 
reduction in the size of the transfer matrix (still in the 
ground state sector) can be judged by comparing the 2nd 
and 3rd columns of Table |l] The above construction has 
enabled us to numerically diagonalize T m for to < 24. It 
is less efficient than the sparse matrix obtained by first 
construction but calculating the nonleading eigenvalues 
is much easier this way. The corresponding correlation 
lengths obtained from the second eigenvalue are given in 
Table IJ 

We have also studied the case when m is not a multiple 
of 3. We recall that in the more familiar case of dimer 
tilings, for even to, there is a one-dimensional height map- 
ping, and accordingly the continuum limit is that of a free 
boson with c = 1. For odd to, however, the height rep- 
resentation has non-periodic boundary conditions, cor- 
responding to a twist operator, which leads to renor- 
malization of the effective central charge to c c s = —2 
|3lj . Returning to the trimer problem, we see that tak- 
ing m (mod 3) ^ also introduces twist operators here, 
and would change the effective value of central charge. 
The entropy per site S m for to (mod 3) = 1 and 2 is 
listed in Table El 



m 


diiTu(r m ) 


dim 2 (T m ) 




£a 


3 


3 


3 


0.37754275 


0.58860147 


6 


21 


6 


0.21764117 


0.95122814 


9 


183 


19 


0.18163298 


1.50340426 


12 


1773 


99 


0.17027036 


2.15126321 


15 


18 303 


672 


0.16557863 


2.83712631 


18 


197157 


5 667 


0.16322214 


3.54463208 


21 


2 189 799 


52 689 


0.16187256 




24 


24 891741 


520 407 


0.16102733 




27 


288 132 303 




0.16046299 





TABLE I: The dimension of the transfer matrices T m , 
dimi(T m ), without symmetrization, and dmi2(T m ), with sym- 
metrization, is shown for different m. Also shown are the 
entropy per site S m and the correlation length £2 defined in 
Eq. (|HSJ - 



VI. NUMERICAL DIAGONALIZATION OF THE 
TRANSFER MATRIX 



T 



3 1 

1 

1 



Two-dimensional isotropic statistical systems with 
(26) short-range interactions at the critical point may exhibit 
invariance under conformal transformations [27|. For a 
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m 


Sm 


m 




4 


0.15682941 


5 


0.13541741 


7 


0.15773925 


8 


0.15031598 


10 


0.15789911 


11 


0.15438343 


13 


0.15806225 


14 


0.15603971 


16 


0.15817988 


17 


0.15687242 


19 


0.15825951 


20 


0.15734824 


22 


0.15831410 


23 


0.15764475 



TABLE II: Entropy per site S m for m(mod3) = 1 and 2. 



cylinder of width to the finite-size corrections to the en- 
tropy per site S m are then of the form [28| 



(27) 



where Soo is the entropy per site in the thermodynamic 
limit to — > oo (i.e., in the infinite plane) and c is the 
central charge (conformal anomaly number) which deter- 
mines the universality class of the problem [27l | . 



(m, m + 3) 




c 


(m, m + 3) 




c 


(m, m + 3) 


Sqc 


c 


(3,6) 


0.16434064 


3.664674 


(4,7) 


0.16728470 


-0.319490 


(5,8) 


0.16121121 


-1.231563 


(6,9) 


0.15282642 


4.456333 


(7,10) 


0.15818038 


-0.041283 


(8,11) 


0.15986633 


-1.167350 


(9,12) 


0.15566128 


4.017785 


(10,13) 


0.15805270 


-0.029334 


(11,14) 


0.15895040 


-1.055394 


(12,15) 


0.15723777 


3.584218 


(13,16) 


0.15829869 


-0.076314 


(14,17) 


0.15871184 


-1.000265 


(15,18) 


0.15786649 


3.314046 


(16,19) 


0.15840838 


-0.111719 


(17,20) 


0.15862738 


-0.968649 


(18,21) 


0.15813526 


3.147735 


(19,22) 


0.15845364 


-0.133847 


(20,23) 


0.15858708 


-0.946406 


(21,24) 


0.15826623 


3.037424 














(24,27) 


0.15833844 


2.957992 















TABLE III: Estimates of entropy Soo and the central charge c from fits using successive pairs of cylinder widths m. 



It is important that estimates of Soo obtained by ex- 
trapolating of data for different values of to (mod 3) have 
to be consistent with each other. Estimates for S^o and 
c obtained by fitting pairs (S m , S m+ ^) to Eq. (|27[1 are 
shown in Tabic ITTll Clearly, there is still some residual 
m-dependence. In general, in conformal field theory, one 
can also have a non-universal 1/m correction to scaling 
term in Eq. (|2"?|) . 



<? —9 



6m 2 



(28) 



with 2 < d < 4. Here we adopt the following strategy: 
we choose trial values of d and S^ and obtain c and 



e from Eq. [2HI by sequential 2-point fits. We find that 
for different to values there is reasonable convergence for 
d = 2.78 ± 0.10 and S^ = 0.158520 ± 0.000015. In table 
IIVI we have listed the values of c and e obtained by two- 
term sequential fits for S m with m = 0, 1 and 2 (mod 3), 
using the fitting form Eq. (|28p. The estimate for central 
charge is c = 2.15 ± 0.2. Our estimate for the effective 
charge is c c s = —0.28 ± 0.02 for to = 1 (mod 3), and 
c off = -0.79 ± 0.02 for to = 2 (mod 3). 



(m, m + 3) 


c 


e 


(m, m + 3) 


c 


e 


(m, m + 3) 


c 


e 


(3,6) 


4.48761407 


-0.8859 


(4,7) 


-0.11255023 


0.0932 


(5,8) 


-0.77423123 


-0.5984 


(6,9) 


2.24822274 


3.8068 


(7,10) 


-0.26174392 


0.4455 


(8,11) 


-0.78830675 


-0.5615 


(9,12) 


1.85242965 


4.9419 


(10,13) 


-0.27727375 


0.4938 


(11,14) 


-0.79454721 


-0.5406 


(12,15) 


1.98106106 


4.4810 


(13,16) 


-0.27266918 


0.4762 


(14,17) 


-0.79175635 


-0.5519 


(15,18) 


2.09438199 


3.9984 


(16,19) 


-0.27314362 


0.4784 


(17,20) 


-0.78906173 


-0.5645 


(18,21) 


2.14525639 


3.7489 


(19,22) 


-0.27959514 


0.5114 


(20,23) 


-0.78903295 


-0.5647 


(21,24) 


2.15734045 


3.6821 














(24,27) 


2.14936018 


3.7310 















TABLE IV: Estimates of c and e from 2-point fits with e/m d correction for d = 2.774 and 5*00 = 0.15852 for m (mod 3) = 0, 1 
and 2. 



To inquire further into the critical behavior of trimer \{m) is some sub- leading eigenvalue (i ^ 1). From con- 
tilings, we can measure as well Sm = log Xi(m), where 
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0.18 
0.17 
0.16 
0.15 




0.02 0.04 



0.06 0.08 

1/m 



0.1 0.12 0.14 



FIG. 8: Convergence of S m with 1/m for m(mod3) = (*), 
m(mod3) = 1 (□), m(mod3) = 2 (o). Continuous lines are 
fits from Eq.EHlthat converges to Soo = 0.15852 for d = 2.774 
and c and e as listed in Table ITVI 



formal field theory we expect to get corrections to follow 
the behavior 29] 



SW = (2A i -|)^ + r,(,» (29, 



6/ m 



where Aj = hi + hi is the scaling dimension of the field 
corresponding to the excited state described by A;. Al- 
ternatively, Eq. (|2~9l may be stated in terms of the corre- 
lation length 



6 



1 



of the excitation on the cylinder. This then reads 

to 

^ = 27a"' 



(30) 



(31) 



i.e., the correlation length is proportional to the cylinder 
width. 

We shall now give numerical evidence that Eq. (|2*9")l 
holds true for a number of excitations (i.e., that the re- 
sulting finite-size estimates for Aj converge as m — > oo). 
This is support for the hypothesis that trimcr tilings are 
conformally invariant. (We note that Kenyon has proved 
[3fj| that dimer tilings are conformally invariant.) 

The first kind of excitation i = 2 in the ground state 
sector (i.e., with the tilings grown on flat base) has al- 
ready been mentioned above. The corresponding correla- 
tion length £ is shown in Tableland plotted against m in 
Fig.0 One finds a slope a = 2 ^ A ~ 0.23, corresponding 
to A 2 ~ 0.69. 

We have also studied the case of 3|m, but with tilings 
grown from a non-planar base. For simplicity we shall 
consider only the two simplest sectors. The first one cor- 
responds to taking the initial state s\ = 1, and Sk = 
for k — 2, . . . ,m, in the notation of Section fV Al The 
associated exponent Ad describes physically the decay 
of the correlation function between a pair of widely sep- 
arated dimer defects in the surrounding soup of trimers. 
The second sector that we shall consider is built from the 



12 



15 



18 



FIG. 9: The correlation length £ plotted against the cylinder 
width m. The dotted line has the slope a = 0.23. 



m 




A D 




Al 


6 


0.11284785 




0.10663883 




9 


0.13675337 


0.511103 


0.13215832 


0.634525 


12 


0.14560151 


0.537242 


0.14225002 


0.684513 


15 


0.14988060 


0.559121 


0.14733680 


0.725974 


18 


0.15229423 


0.577165 


0.15029069 


0.761116 


21 


0.15379778 


0.592680 


0.15217226 


0.791808 



TABLE V: Entropy per site for tiling of non-planar bases and 
estimates of scaling dimensions with l/m d correction with d — 
2.774 correction. 



initial state s\ — 1, S2 — 2, and Sk — for k = 3, . . . , to. 
The associated exponent Al corresponds to a pair of L- 
shaped defects (or closely bound compounds of a dimer 
and a monomer). 

The estimates for Ad and Al, obtained by adding a 
non- universal l/m d correction to Eq. I|29|) . are shown in 
Table [V] The dimer defect is the most relevant. In the 
continuum limit, described by the two-dimensional height 
field, the exponents Ad and Al would describe the de- 
cay of vortex-vortex correlation functions. Their final 
values appear to be of the order Ad = 0.62 ± 0.1 and 
Al = 0.82 ±0.1, but clearly the words of caution on the 
slow convergence rate made above when discussing the 
extrapolations of c are also applicable here. 



VII. CORRELATION FUNCTIONS 

Given a trimer tiling, we assign to each lattice face a 
state s = 1,2, ... ,6 according to its position in the tiling. 
Our convention is that horizontal trimers are labeled as 



1 2 3 



and vertical trimers as 



Then we can define the correlation function Gjj (x, Y) 
as the probability that the face at X is in the state s = i 
and the face at Y is in the state s = j (with 1 < i,j < 
6). In the thermodynamic limit, translational invariance 
implies that Gj^^X, Y) is only a function of (X — Y), 
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and we write it as Gy (R) with R = X Y. However, 
all the correlation functions are not independent. For 
instance, we can express the function G\j (R) in terms of 
the function G2j(R) as 



G lj (R) = G 2j (R + e x ). 



(32) 



Using such equations, one can express the functions where 
i or j take values 1,3,4 or 6 in terms of the functions 
G22,G25 and G55. Again G22 or G55 are related by the 
symmetry in the horizontal and vertical direction. This 
only leaves the functions G22 and G25. The conditional 
probability that the site R is occupied horizontally, given 
that the origin is in state 2, is written 1/2 + /(R), so 
that the probability of finding a vertical trimer at R is 
1/2 - /(R). We therefore have 

G 22 {TL)+G 2 2(R + e x ) + G 22 (R-e x ) = i+/(R)(33) 

G 25 (R) + G 25 (R + e y ) + G 25 (R -ey) = ±- /(R)(34) 

The above equations can be solved as linear inhomoge- 
neous equation in G22 and G25, and can be solved if /(R) 
is known for all R. Thus we have only one unknown func- 
tion /(R). Note that Eq. (|33|l is a difference equation in 
the x-coordinate (the y-coordinate does not change). The 
general solution reads 

G 22 (X,Y) = G° 22 (X) + A Y sin k Q X + B Y cos k X (35) 

where fco = 27r/3. But A and B must be zero, since the 
function G22 (X, Y) is known to decay to a constant value 
1 /6 for large X. Thus we can determine G22 and similarly 
G25 uniquely in terms of /(R). 



VIII. MONTE-CARLO SIMULATIONS 

We have already introduced the local move of a block 
of three vertical trimers replacing a block of three hor- 
izontal trimers, or vice versa (see Fig. 0J. It has been 
shown that on a square lattice the above move is ergodic, 
i.e. , all possible configurations can be reached from any 
initial configuration |22( . Starting from any initial config- 
uration, say the standard configuration with all trimers 
vertical, repeating the local operation at randomly chosen 
sites generates a random trimer tiling. 

In our Monte Carlo simulation, we start with a L x 
L lattice fully packed with all trimers vertical (L = 
15,45,60,90). In one Monte Carlo step, we randomly 
select one of the trimers, and check if can be flipped us- 
ing the move of Fig. If yes, the flip is made, if no, the 
move is rejected, and another site is selected. We discard 
the initial 10 5 steps. Once the steady state is reached, we 
calculate the different correlation functions in the steady 
state. We generated data for over 10 5 Monte-Carlo steps 
per site (MCS). 

To verify that our sampling produces unbiased results, 
we have used it to compute the average height on each 
of the nine different sublattices. These can also be com- 
puted analytically as follows. First, we define the ab- 
solute values of the heights by setting the mean height 



(h)i = on sublattice 1. To compute (h)2 on sublat- 
tice 2, we first note that the height difference between 
site A, belonging to sublattice 1, and site B, belonging 
to sublattice 2, (i.e., along the horizontal edge joining 1 
and 2 in Fig. Ufa)) can either be 1 if the edge overlaps 
with the boundary of a horizontal or vertical trimer, or 
it can be —2 or 2ui — uo 2 if it is one of the internal edges 
of a vertical trimer. Hence the mean value of — is 
(l/6)(4+2w — to 2 — 2) = ~uo 2 /2, and since B was arbitrary 
(h)2 = —uj 2 /2. In a similar way, we find that 



h 


1 = 


(h) 


5 = 


(h 


9 


= 


h 


2 = 


(h) 


6 = 


(h 


7 




h 


3 = 


(h) 


4 = 


<h 


8 


= W 



We have found these values of (h)2 and (h)3 to be in 
excellent agreement with the numerical results. 
To compute the height-height correlation function 



H(x,y) = (\h(X 



x,Y + y)-h(X,Y)\' 



(36) 



we averaged the data over all positions (X,Y). The re- 
sults are shown in Fig. 1101 For next nearest neighbors, 
height difference squared is 1 with probability 2/3, 4 and 
7 each with probability 1/6. Hence if (1,0) = H(0, 1) = 
5/2 which has been verified against the results of the sim- 
ulations. For system size L, H(x,y) varies for large r as 



H(x, y) « A 2 log [{BLf {sin 2 (^) + sin 2 (^) } + 1 



(37) 



where A 2 = 0.69 and B = 4. Since H(l, 1) ^ H(l, -1), 
there should also be a term depending on sin 2 (7ra/L) — 
sin 2 (iry / L) , but this should decay faster for large r, 
and we have not included it in our fits. For large 
r = y/ 1 x 2 + y 2 it is easy to see that 



H(r) m 2A 2 log(r) 
or equivalently we have 

c - G(r) w A 2 log(r) 



(38) 



(39) 



where c is a constant and C(r) = (h(r)h(O)}. From our 
estimates of A2 = 0.69 and a = 0.23, from Eq. I|31|) it is 
seen that A2 = l/(27ra) in agreement with the theory of 
conformal invariance. 

In Fig.^l we have plotted f(R, 8) for three different di- 
rections 9 = 0, 7r/4, it/2. We used a lattice of size 90 x 90, 
and 10 5 MCS to get the data. In each case, we see that in 
each direction, the function f(R, 9) decreases as a power, 
/ ~ R~ x , with x w 1.5. However it is observed that the 
effective exponent decreases with R and it is difficult to 
make meaningful estimates due to large fluctuations in 
the data. 



IX. CONCLUSIONS 

We have studied the problem of tiling the plane with 
trimers, and seen that it differs from the well-known 
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(dimers and L-shaped trimers). 




FIG. 10: Correlation function H(x,y) for L=60. 
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FIG. 11: \f(R,6)\ as functions of R for 9 = (•), = tt/4 
(o), = 7r/2 (□). Straight lines have slope -1.5. 



dimer tiling problem in a number of ways. First, 
the number of sectors in the transfer matrix grows 
exponentially rather than linearly with system size. 
Second, the natural height mapping has a target space 
dimension which is two rather than one. If the tiling 
problem is conformally invariant, one could therefore 
reasonably expect it to have central charge c = 2. To 
check this hypothesis, and the inference on c, we have 
performed extensive numerical transfer matrix calcula- 
tions, on cylinders of widths m < 27 lattice spacings. 
The finite-size corrections exhibit unexpectedly large 
correction-to-scaling terms, making difficult the precise 
estimations of the critical exponents. The results are 
however clearly in favor of conformal invariance, and this 
is confirmed by the logarithmic form of the height-height 
correlations as obtained by Monte-Carlo simulations. 
Our estimate c = 2.15 ± 0.2 is in marginal agreement 
with the expectation c = 2, but we cannot definitely 
rule out a more complicated behavior. We have also 
measured numerically a number of other exponents, in 
particular the leading thermal scaling dimension and the 
scaling dimensions related to pairs of geometrical defects 
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APPENDIX A: IRREDUCIBLE WORDS FOR A 
GENERAL TILING PROBLEM 

The working of section IV can be adapted to a more 
general tiling problem in which there are two types of 
straight tiles: horizontal p-mers (i.e., of size px 1 ele- 
mentary squares) and vertical g-mers (of size 1 x q). We 
suppose for the rest of this appendix that p > 3 and that 
q > 2. 

The words describing the sector decomposition of the 
transfer matrix are then made of q different letters. An 
irreducible word is one in which there is no substring of 
p consecutive equal letters. For 1 < k < p — lwe then 
define Nk (n) as the number of irreducible words of length 
n in which the last k letters are all equal, but in which 
the last k + 1 letters are not all equal. 

The initial conditions for the recursion relations are by 
convention 



N k (p -l) = lforl<fc<p-l 



(Al) 



This definition of Nk(n) does not yet take into account 
the multiplicity due to the q-dependent number of ways 
that one may chose the last p — 1 letters. Denoting this 
multiplicity M&(p — 1) one finds 



M k (p-l) 
M p _i(p-1) 



p—i—k 



(q - 1) for 1 < k < p - 2 



(A2) 



Indeed, one needs to choose the first p — 1 — k letters ar- 
bitrarily and then complete the word with a single letter, 
chosen different from the last one chosen arbitrarily. Note 
that the sum of all multiplicities is Ylk=i Mk(p — 1) = 



as expected. 



For any n > p — 1 one then has the recursion relations 



p-2 

JV x (n+l) = («-l)Xy-*JVfc(n) + g 3 - J 'JVi 
k=i 

N k (n + 1) = qN k -i(n) for 2 < k < p - 2 



(n) 



A^ p _i(n + 1) = (<f-l)JVp_ 2 (n) 



(A3) 



This generalizes Eq. (JHJ ■ Equivalently, the recursion rela- 
tions can be written in matrix form by defining the vector 
N(n) with elements Nk(n). One then has N(n + 1) = 
TN(n), where T = {TV} is a matrix with elements that 
can be read off from 1 A3(l . Iterating this, one finds 



N(n) 



-Pl = ST> n+1 - p S~ 1 T, 



(A4) 



where 1 = [1,1,..., 1]* and S is the matrix of eigen- 
vectors of T that turns T into diagonal form through 
D = S _1 TS. Writing this out for given values of p and 
q will yield an explicit solution generalizing Eq. (J7J. 



The total number t n of irreducible words of length n is 
then obtained by taking into account the multiplicities: 



p-i 

t n = ^M k {p-l)N k (n) 

k=l 



(A5) 



This generalizes Eq. ©. 

To obtain the total number of words that are irre- 
ducible to the empty word, we first define Gg as the num- 
ber of irreducible unfactorizable words with initial letter 
£ = 1, . . . , q. Introducing also 



n=0 \i=2 , 

the generalization of Eq. HUH reads 
G t = t{G Q l) p - X 



(A6) 



(A7) 
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The corresponding generating function g(x) then satisfies 



g(x) [l-{q-l)g(x)r l =^ 



(A8) 



The generating function for all words that are reducible 
to the empty word is then 



(A9) 



and reads explicitly 

H{x) 



1 



1 - qg(x) 



(A10) 
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